Breathing dynamics in heteropolymer DNA 
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While the statistical mechanical description of DNA has a long tradition, renewed interest in 
DNA melting from a physics perspective is nourished by measurements of the fluctuation dynamics 
of local denaturation bubbles by single molecule spectroscopy. The dynamical opening of DNA 
bubbles (DNA breathing) is supposedly crucial for biological functioning during, for instance, tran- 
scription initiation and DNA's interaction with selectively single-stranded DNA binding proteins. 
Motivated by this, we consider the bubble breathing dynamics in a heteropolymer DNA based on a 
(2-f l)-variable master equation and complementary stochastic Gillespie simulations, providing the 
bubble size and the position of the bubble along the sequence as a function of time. We utilize new 
experimental data that independently obtain stacking and hydrogen bonding contributions to DNA 
stability. We calculate the spectrum of relaxation times and the experimentally measurable autocor- 
relation function of a fluorophore-quencher tagged base-pair, and demonstrate good agreement with 
fluorescence correlation experiments. A significant dependence of opening probability and waiting 
time between bubble events on the local DNA sequence is revealed and quantified for a promoter 
sequence of the T7 phage. The strong dependence on sequence, temperature and salt concentration 
for the breathing dynamics of DNA found here points at a good potential for nanosensing applica- 
tions by utilizing short fluorophore-quencher dressed DNA constructs. 

Key words: DNA denaturation; biomolecules; fluorescence correlation spectroscopy; master equa- 
tion; stochastic simulation 
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I. INTRODUCTION 



Textbook pictures of double-stranded DNA molecules 
may lead one to believe that the Watson-Crick double- 
helix represents a static geometry. As a matter of fact, 
even at room temperature DNA opens up intermittent 
flexible single-stranded domains, so-called DNA-bubbles. 
Their size typically ranges from a few broken base-pairs 
(bps), increasing to some 200 broken bps closer to the 
melting temperature P, [3, i, 3 . The stability of DNA 
is characterized by the two free energies ehb for breaking 
the Watson-Crick hydrogen bonds between complemen- 
tary AT and GC bps, and the ten independent stack- 
ing free energies est for disrupting the interactions be- 
tween neighboring bps; at 100 mM NaCl concentration 
and temperature 37 °C it was found that ehb = l.OfcsT 
for a single AT and 0.2kBT for a single GC-bond (at 
T = 37°C, ksT = 0.62kcal/mol). Under the same con- 
ditions the weakest (strongest) stacking interaction was 
found to be the TA/AT (GC/CG) with free energies 
Est = -O.OfcsT (-4.1fcBT) H]. In addition, the initiation 
of a bubble in an unperturbed DNA molecule, creating 
two interfaces between single-stranded bubble and vici- 
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nal double- helix at the zipper forks, is associated with 
an activation factor (Tq ~ 10"^ . . . 10~^ [ijSSHi related 
to the ring-factor ^ used in Q and below. That is, 
despite the rather low free energy for breaking the bps, 
the high bubble nuclcation barrier guarantees that below 
Tm bubbles are rare and well separated, in particular, 
under physiological conditions. However, once a bubble 
opens, since typical free energies are of the order ksT lo- 
calized denaturation zones can open up, predominantly 
in AT-rich regions [l|,|3,ll. These DNA-bubbles fluctu- 
ate in size, the DNA-breathing. It has been demonstrated 
recently by fluorescence correlation methods that DNA- 
breathing can be probed on the single molecule level, 
revealing a multistate kinetics of stepwise (un)zipping of 
bps with a bubble lifetime ranging up to a few millisec- 
onds 0. 

Theoretically, based on the statistical mechanical 
Poland-Scheraga model DNA-breathing has been de- 
scribed in homopolymer DNA in terms of a continuous 
Fokkcr-Planck equation [l^, and through a stochastic 
Gillespie scheme [ll|. A discrete master equation ap- 
proach was developed in [l^, , including the coupled 
(un)binding dynamics of selectively single-stranded DNA 
binding proteins. Continuous and discrete approaches 
are compared and studied in [l3|. Heteropolymer DNA- 
breathing was considered in a reduced one-variable ap- 
proach using a random energy model (isj . 

Here, we develop a full (2 + l)-variable approach to 
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breathing in heteropolymer DNA that allows us to study 
the sequence dependence of the dynamics, through the 
initiation and the stochastic motion of the two forks of 
an open DNA bubble. Two approaches are used: the 
stochastic motion is obtained by generating stochastic 
(Gillespie) time series from which equilibrium distribu- 
tion as well as autocorrelation functions arc obtained. 
We also use the corresponding master equation in or- 
der to calculate the complementary ensemble-averages; 
excellent agreement is found between time-averages and 
ensemble averages. Novelties in our study include: (i) 
We study the full dynamics (2+1 variable problem) of a 
heteropolymer region of arbitrary (not just random) se- 
quence; (ii) we compare our results to the fluorescence 
correlation spectroscopy (FCS) experiments in Q using 
the directly measured DNA parameters in 0] (see be- 
low); (iii) Recently, for the first time, the (two) hydrogen 
bond energies, and (ten) stacking interactions character- 
izing DNA stability within the Poland-Scheraga model 
were separately determined d, |1] ; these stability param- 
eters are utilized in our study. Among the consequences 
of these new results compared to previously used param- 
eters are the more pronounced sequence dependence 
and the fact that Watson-Crick and stacking interactions 
can completely be separated as required when studying 
internal bubble dynamics (a bubble involving m broken 
Watson-Crick bonds and m + I broken stacking interac- 
tions). 



Based on this new approach, wc study the question 
of transcription initiation at the TATA motif of the bi- 
ological sequence in the bacteriophage T7 promoter se- 
quence. Using the newly obtained stacking parameters 
from [§| , we demonstrate the delicate dependence of both 
the equilibrium opening probability as well as the breath- 
ing dynamics on the sequence dependence of the stacking. 
While in our model the opening times of bubbles only 
marginally depend on their position along the sequence, 
the recurrence frequence of bubble events is much more 
sensitive to the position. The latter might therefore be a 
clue toward the understanding of transcription initiation. 



This paper is organized as follows: In Section |TT] we 
describe the DNA bubble dynamics in terms of the rel- 
evant transfer coefficients. In Section IIII Al a stochastic 
scheme based on these transfer coefficients in terms of 
the Gillespie algorithm is introduced. In Section IIIIBI a 
complementary master equation scheme is described. In 
Section [IVI we apply our two complementary formalisms 
to (i) the experimental constructs in Ref. [9j ; (ii) the T7 
phage promoter sequence, (iii) We show a strong depen- 
dence on sequence, temperature and salt concentration 
and demonstrate the good potential for nanosensing ap- 
plications. Technical details necessary for the introduc- 
tion of our model appear in a separate publication p^ . 



fluorophore 




quencher 



FIG. 1: Clamped DNA domain with internal bps x = 1 to M, 
and tag position xt- The DNA sequence enters through the 
statistical weights list (a;) and uwa{x) for disrupting stacking 
and hydrogen bonds respectively. The bubble breathing pro- 
cess consists of the initiation of a bubble and the subsequent 
motion of the forks at positions xl and xr. 



II. GENERAL MODEL AND TRANSFER 
RATES 

With typical experimental setups Q in mind, we con- 
sider a segment of double-stranded DNA with M inter- 
nal bps, that are clamped at both ends, i.e., the bps 
X = and x = M + 1 are always closed (Fig. [1]). The 
heteropolymer character of the problem enters via the 
position-dependence of the statistical weights Mhb(a;) = 
exp{ehb(-T)/(fcBr)} for breaking the hydrogen-bonds of 
the bp at position x, and Ust{x) = cxp{est(a;)/(fc_BT)} for 
disrupting the stacking interactions between bps x—1 and 
X] est (a;) and e\±,(x) are the corresponding free energies, 
which in general have energetic as well as entropic contri- 
butions. Due to the high free energy barrier for bubble 
initiation (^ <IC 1, see below), opening and merging of 
multiple bubbles are rare events, such that a one-bubble 
description is appropriate [l^. The positions xl and xr 
of the zipper forks respectively correspond to the right- 
and leftmost closed bp of the bubble; these are stochastic 
variables whose time evolution characterizes the bubble 
dynamics. Note that writing the Boltzmann factors for 
the free energies as exp {AG/(fcBT)}, a positive AG de- 
notes an unstable bond. 

In terms of xl and bubble size m = xr — xl — 1, the 
bubble partition factor is (m > 1) 

^/ XL+rn XL+m+1 

x—xl+I x—xl+I 

completed by 3f{m = 0) = 1. Here, ^' = 2'^^, where £^ 
10"^ is the ring factor for bubble initiation from Ref. l8| . 
For a homopolymer £ is related to the cooperativity pa- 
rameter cTo « 10"^ 0, 01 by (To = £,exp{est/{kBT)} [1]. 
For the entropy loss on forming a closed polymer loop 
we assign the factor (1 -I- m)~^ |y, [13] and take c = 1.76 
for the critical exponent [l3|. Note that a bubble with 
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m open bps needs breaking of m hydrogen bonds and 
m + 1 stacking interactions, see Eq. ([1]). The equiUbrium 
probabihty for finding a bubble with a given xl and m 
is 



(2) 



Below we impose the detailed balance condition when 
introducing the rates to guarantee that Pcq{x]^,m) is in- 
deed reached for long times. 

Let us proceed by introducing the transfer (rate) coef- 
ficients for the bubble dynamics. For the left zipper fork 
we define {xl , m) as the transfer coefficient for the pro- 
cess XL —> XL -I- 1, corresponding to bubble size decrease, 
and t~[{xLi m) as the transfer coefficient for xl — > xl — 1 
(bubble size increase). For the right zipper fork we sim- 
ilarly introduce ^^{xl, m) for a;^ — > -|- 1 (bubble size 
increase) and t^ {xL, m) for xji — ^ xr — 1 (bubble size de- 
crease). In addition for the transition m = ^ m = 1, 
i.e., for the initial bubble opening process occurring at 
position XL, we introduce ^^(xl), and for the bubble 
closing process TO=I^m = Owe employ ^^(xl). Note 
that t(ij(xL) and ^^(xl) correspond to closing or opening 
of the bubble at position x = xl + f . Due to the clamp- 
ing we require that xl > and xji < A/ + I, and we 
therefore introduce reflecting conditions 

t^ixL ^ 0,rn) = t^{xL,m ^ M - Xl) ^ (3) 

(also, tJ(xL = —1, m) = and t^(xL, m = M— xl + 1) = 
for 7Ti = 2, i\/ + 1 for completeness). 

Let us consider explicit forms for the transfer coeffi- 
cients. For bubble size decrease we take 

tJ(xL,m)|„i>2 = t~^{xL,m)\ra>2 = /C(to)/2 (4) 

for the left fork and right forks, respectively. We above 
defined the m-dependent rate coefficient 



A^(to) = km ^ . 



(5) 



As in previous studies, this expression imposes the hook 
exponent /x, related to the fact that during the zipping 
process not only the bp at the zipper fork is moved, but 
also part of the vicinal single-strand is dragged or pushed 
along [l^, [l^, [l^l . One would expect that the hook expo- 
nent is only relevant for larger bubbles, and we put /i = 
in the remainder of this work, mainly focusing on T well 
below where the bubbles sizes are small. The rate k 
characterizes a single bp zipping. Its independence of x 
corresponds to the view that bp closure requires the dif- 
fusional encounter of the two bases and subsequent bond 
formation; as sterically AT and GC bps are very simi- 
lar, k should not significantly vary with bp stacking, k 
is the only adjustable parameter of our model, and has 
to be determined from experiment or future MD simula- 
tions. The factor 1/2 is introduced for consistency with 
previous approaches (T^ . IT3| . We note that, in princi- 
ple, an x-dcpendcnce of k can easily be introduced in our 



approach by choosing different powers of the statistical 
weights entering the rate coefficients such that they still 
fulfill detailed balance. 

Bubble size increase is controlled by 

t£(xL,m) ^ /C(m -M)wst(xL)uhb(xL)s(m)/2, 
tJ(xL, m) = /C(m + l)wst(xH + l)uhb(a;_R)s(m)/2,(6) 

for m > 1, where 



3(m) = {(l + m)/(2 + 7n)}^ 



(7) 



For m > 1 we thus take the rate coefficients for bubble 
increase proportional to the Arrhcnius factor Ust^hb = 
exp{(ehb + est)/[fc_B2^]} multiplied by the loop correction 
s(m). Note that an unzipping event on average involves 
the motion of one more open base-pair compared to a zip- 
ping event, and the transfer coefficients above are there- 
fore proportional to K.{m + 1). Finally, bubble initia- 
tion and annihilation from and to the zero-bubble ground 
state, TO = <-> 1 occur with rates 

^g^^l) = k^'s{0)ust{xL + l)uhh{xL + l)UstixL + 2) 



tp(xL) 



(8) 



with the bubble initiation factor ^' included in the ex- 
pression for tg. Note that t^, in contrast to the opening 
rates for to > 1, is proportional to an Arrhcnius- factor 
involving two units of stacking free energy. The annihila- 
tion rate Iq^xl) is twice the zipping rate of a single fork, 
since the last open bp can close either from the left or 
right. The rates t together with the boundary conditions 
fully determine the bubble dynamics. 

We see that the rates t^, t^, and t^ are chosen such 
that they fulfill the detailed balance conditions: 



^li^L - 1,™+ l)Pc(iixL - 1,toH 

= tl{xL,m)Pcq{xL,m), 
tJi{xL,m + l)PcqixL,m + 1) 

= t+{xL,m)Pcq{xL,rn), 

t(i;(xL)Poq(0) t^{xL)Pcq{xL, 1) 



1) 



(9) 



These conditions guarantee relaxation towards the equi- 
librium distribution P°'^{xl, m), see Eq. In the next 
two sections we use the above explicit expressions for 
the transfer coefficients and describe the DNA breath- 
ing dynamics pursuing two complementary approaches: 
the stochastic Gillespie scheme (Section IIII Ap and the 
master equation f Section IIII Bp . 



III. DYNAMIC APPROACHES TO 
DNA-BREATHING 

A. Gillespie approach 

In this section we use the Gillespie algorithm together 
with the explicit expressions for the transfer coefficients 
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introduced in the previous section in order to generate 
sequence specific stochastic time series of breathing bub- 
bles. In particular we show how the motion of a tagged 
bp is obtained. 

To denote a bubble state of m broken bps at position 
XL we define the occupation number b[xL,'m) with the 
properties b{xL, m) = 1 if the particular state {xl, ti} is 
occupied and b{xL, tn) — for unoccupied states. For the 
completely zipped state m = there is no dependence on 
Xl, and we introduce the occupation number 6(0). The 
stochastic DNA breathing then corresponds to the near- 
est neighbor jump processes in the lattice of permitted 
states [16|. In the Gillespie scheme, each jump away from 
the state {xl, m} (i.e., from the state with b{xL,'m) ~ 1) 
occurs at a random time r, and in a random direction 
to one of the nearest neighbor states. This stochastic 
process is governed by the reaction probability density 
function [nl, [m, m 

P(T,/i,!/) = t;;(a;i,m)exp I -T^X.'^,{xL,m) j . (10) 

More explicitly, for a given state {xl, m) the joint proba- 
bility density (fTO]) defines after what waiting time t after 
the previous random step the next step occurs, and in 
which reaction pathway, ly G {G,L,R}, fj, £ {+/—}■ In 
the present case, v and /i denote x-dependent zipping 
or unzipping of a bp at the left of right zipper fork. A 
simulation run produces a time scries of occupied states 
{xL,'m} and how long time r = Tj {j ~ I, ...,iV, where 
N is the number of steps in the simulation) this partic- 
ular state is occupied. This waiting time r, in particu- 
lar, according to Eq. PU)) follows a Poisson distribution 
[23j . Note that the waiting times governed by (fTO|) vary 
widely, as the reaction rates occur in the exponential (in 
particular, the bubble initiation with the ^-factor has a 
long characteristic time scale). The fact that the Gille- 
spie scheme uses the weighted reaction time scale instead 
of fixed simulation time steps makes this algorithm very 
efficient. 



takes only values or 1) over region Ml, i.e., 

Ht)^ b{xL,m). (12) 

From the time series for I{t) one can, for instance, calcu- 
late the waiting time distribution ^/'(r) of times spent in 
the / = state, as well as the survival time distribution 
(/)(t) of times in the 1=1 state. Explicit examples for 
iP{t) and (/)(t) are shown in Section [TVl 

The probability that the tagged bp is open becomes 

1 ^ 

PGiti) = J^T.^lHtj), (13) 

where tj = ^'^^ long times the explicit con- 

struction of the Gillespie scheme together with the de- 
tailed balance conditions guarantee that Paitj) tends to 
the equilibrium probability, i.e., that Paitj — > c>o) = 
J2xL mmi^'^'^i^L^m), where P^i(a;L,m) is given in 
Eq. @. 

2. Tagged base-pair autocorrelation function 

The autocorrelation function for a tagged bp is ob- 
tained through 

At{xT,t) = mm-(mf 

= I{t + t')I{t')dt' 1^ I{t')dt'^ (14) 

which for long sampling times T converges to the en- 
semble average, Eq. (jlSp . from the master equation (in- 
troduced in the next section). The function At{xT,t) 
corresponds to the quantity obtained in the fluorescence 
correlation experiment of Ref . [§] . 

B. Master equation formulation 



1. Tagged bp survival and waiting time densities 

Motivated by the experimental setup in we study 
the motion of a tagged bp at a; = xt , see Fig. [TJ In the 
fluorescence correlation experiment fluorescence occurs if 
the bps in a A-neighborhood of the fluorophore position 
XT are open Measured fluorescence time series thus 
correspond to the stochastic variable I{t), with the prop- 
erties I{t) = 1 if at least all bps in [xt — A, xt + A] are 
open, and I(t) = otherwise. Thus, if / = 1 we are in 
the phase space region defined by 

Ml : {0 ^ XL ^ XT-A-1, XT-XL + A < to < M-xl}- 

(11) 

Conversely, / = corresponds to the complement RO of 
Ml. The stochastic variable I{t) is then obtained by sum- 
ming the Gillespie occupation number b{xL, to) {b(xL, to) 



Complementary to the stochastic simulations of DNA- 
breathing detailed in the preceding section we here intro- 
duce a master equation for the joint probability density 
P{t) = P{xL,m,t;x']^,m' ,0) that at time t the system 
is in state {xL,m} and that it was in state {x^,to'} at 
t = 0. The master equation, which is equivalent (in the 
sense that it produces the same averaged quantities) to 
the Gillespie scheme, can be formally written as 

|p(t) = WPit), (15) 

where the explicit form of the matrix W is given in 
terms of the rate coefficients from the previous section 
in Ref. [lB|. A standard approach to the master equa- 
tion is the spectral decomposition [l^, l25l | 

Pit) = Y,CpQpeM-Vpt)- (16) 
p 
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The coefficients Cp are obtained from tfie initial condition. 
Inserting Eq. pB]) into Eq. ([TS]) produces tlie eigenvalue 
equation 

l^Qp = -VpQp- (17) 

From the eigenvalues rjp and eigenvectors Qp of Eqs. (fT7)) . 
any quantity of interest can be constructed. 

C. Dynamic quantities for a tagged bp 

The waiting time density ij.'it) and the survival time 
density 4>{t), as obtained in a Gillespie scheme, corre- 
spond to the first passage problem to start from an initial 
state with 7=1 (7 = 0) and passing to / = (/ = 1). 
It is discussed in detail in Ref. [lB| how these quantities 
can be obtained from the master Eq. ([TS]). 

The equilibrium autocorrelation function 

A{xT,t)^{i{t)m)-i{i)r (18) 

is a measure for the relaxation dynamics of the tagged 
bp. This can be seen from the identity 

1 1 

(/(t)/(0)) = ^ ^ Ip{I, t; I', 0)1' = p(l, t; 1, 0), (19) 

7=0 /'=0 

where p(l,t; 1,0) is the survival probability density that 
I{t) = 1 and that /(O) = 1 initially. Using the fact that 
p(l, t; 1, 0) is obtained by summing P{xl, m, t; x'j^,m', 0) 
exclusively over region Ml we obtain 

p(l,t;l,0)= J2 P{xL,m,t;x'L,m\0). (20) 

Combining this result with Eq. ([Tn)) . the spectral decom- 
position (|16p . and assuming that we initially are at equi- 
librium: P(xL,TO,0;x'^,m',0) = Sr,im'SxLx'j-PcqixL,m), 

the autocorrelation function p8|) can be rewritten as 

AixT, t)=J2 [Tp{xT)f eM-t/rp), (21) 

with relaxation times Tp = l/rjp, and where 

rCT-A-l M-xl 

Tp{xT)= J2 Qp(^L,m). (22) 

xj^—O rn—xj-—XL-\-A 

For long times, i.e., when the time average is long enough, 
A{xt, t) agrees with At{xT, t) given in Eq. ((M]) as will be 
illustrated in the next section. We can rewrite the corre- 
lation function according to the spectral decomposition 

A{xT,t)= J dTCXp{~t/T)fixT,T), (23) 

where we introduced the weighted spectral density 

fixT,T)=Y,[Tp{xT)?SiT-Tp). (24) 

p#0 



This relaxation time spectrum directly provides the spec- 
tral content of the relaxation behavior of the DNA- 
bubble, and sometimes a better (but equivalent) visual- 
ization of the system than the autocorrelation function. 



IV. RESULTS 

In this section we apply our two complementary for- 
malisms to study the behavior of (i) the designed DNA 
constructs used in the experiments of Ref. [9|, and (ii) 
the T7 phage promoter sequence. 



A. Comparison to experimental results 

In Fig. [2] the autocorrelation functions At{xT,t) for 
the sequence AT9 from Q are shown for various tem- 
peratures T. The data were scaled by k such that the 
curves coincide where A{t) = 1/2. The strong scatter at 
short times is mainly ascribed to quantum transitions in 
the fluorophore 0, [2^ . The lower graph shows the tem- 
perature dependence of the characteristic zipping time, 
1/k. Individual autocorrelations for three temperatures 
are compared in Fig. [3l 

In the combined autocorrelation plot. Fig. [21 the black 
line shows the predicted behavior of A{xt, t), calculated 
by numerical solution of the eigenvalue Eq. (|17[) by help 
of Eq. dsn). Stability parameters from 0] for T = 49°C 
and 100 mM NaCl concentration were used. As in the 
experiment we assumed that fluorophore and quencher 
attach to bps xt = 17 and xt + 1, and that both are 
required open to produce a fluorescence signal (the out- 
ermost GC-pairs in the sequence given in figure [2] were 
taken as clamped, i.e., labeled as a; = and x = M + 1). 
From the scaling plot, we calibrate the zipping rate as 
k = 7.1 X 10*/s for T = 49° in good agreement with 
the findings from Ref. @. The calculated behavior re- 
produces the data within the error bars. The green 
curve corresponds to the ME result for T = 33° C, show- 
ing more pronounced deviations from experimental data. 
Notice that for lower temperatures the relaxation time 
distribution /{xt^t) becomes narrower (Fig. [2] inset). 
Thus, our model predicts that the dynamics for smaller 
temperatures involve fewer modes, which is in contrast 
to the experimental data that have a broad, multimodal 
behavior also for low temperatures. 

The individual behavior of the autocorrelation is dis- 
sected in Fig. [3] for three temperatures spanning the full 
T-range probed in the fluorescence experiments. Note 
the good quality of the match between experimental data 
and model prediction for the highest temperature (49° C). 
This temperature is already close to the denaturation 
temperature of the bubble domain of the AT9 construct 
(the contribution of the longest relaxation time in the 
rather broad spectrum of relaxation times is consider- 
ably larger than the three previous ones). The tendency 
of overestimation of the slope in the autocorrelation func- 
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FIG. 2: Top: Autocorrelation function At{xT,t) at various 
temperatures T measured for the sequence AT9 from 0] at 
100 mM NaCl. The sequence is indicated in the figure, where 
the four lower case t's symbolize a small bulge loop. The 
full lines show the results from the master equation based on 
the DNA parameters from Q. Inset: Relaxation time spec- 
trum /(r) = f{xT,T), showing broadening with increasing 
temperature. Bottom: Characteristic zipping time 1/k as a 
function of temperature in an Arrhenius plot. The line shows 
a least squares fit to an Arrhenius law r oc exp(j4/r) with 
A = 2.6x lO^K. 



tion by our model at lower temperatures is obvious for 
curves at 22°C and 33 °C, while the experimental slope 
remains almost constant over this T-range. 

We expect three effects to contribute to the deviations 
by broadening the relaxation time spectrum, i.e., lower- 
ing the free energy of the system: 

(i) In the present fluorescence correlation spectroscopy 
experiments, two contributions superimpose to produce 
the fluorescence signal [2^ : The diffusional motion of the 
molecule carrying the fluorophorc in and out of the con- 
focal volume, and the actual breathing dynamics. With- 
out the breathing the autocorrelation function takes the 
form A{t) ~ 1/(1 -t- t/Tr)) (for a narrow beam waist), 
where td = w'^ /{AD) « 150ms, with w being the linear 
size of the beam waist and D is the diffusion constant of 




0.1 1 10 100 



t [1/k] 

FIG. 3: Individual autocorrelation functions At{xT,t) for 
three different temperatures (22° C, 33° C, 49° C) spanning the 
T-range probed in the fluorescence experiment. While for the 
highest temperature, the match between data and theory is 
very good, deviations occur at lower temperatures. Reasons 
for these deviations are discussed in the text. Note that the 
data for 33° C are shifted vertically by 0.2, the ones for 22° C 
by 0.4. 



the construct. In [9j the pure diffusive contribution was 
eliminated by performing a separate experiments with 
the quencher being removed (measuring the solely diffu- 
sive contribution), and dividing out this result from the 
signal. However, as the quencher is removed the diffusion 
constant of the construct is slightly changed. In order to 
roughly account for this fact, the blue curve shown in 
Fig. [2] was obtained by a 3% reduction of the diffusion 
time T£). Note that the agreement of the blue line with 
the data is excellent. This underlines the sensitivity of 
the DNA-breathing single molecule data, pointing toward 
potential dynamic methods to calibrate both k and AG. 

(ii) It is very likely that the presence of the fluorophore 
and quencher molecules destabilizes the DNA — despite 
the short stalk through which fluorophore and quencher 
are attached — by altering the Watson-Crick and stacking 
interactions. The resulting decrease of the stacking free 
energy therefore is expected to effect a lower free stacking 
energy in comparison to the undressed DNA, for which 
the stability data are measured and which are used in 
our model. 

(iii) Finally, our present model does not take into ac- 
count the entropic contributions due to the degrees of 
freedom of the ffuorophore/quencher pair, i.e., the fact 
that for bigger bubbles the ffuorophore/quencher pair 
has more freedom to diffuse around and rotate. To ap- 
proximately account for this we would change the par- 
tition from 3f{xL,m) to riFQ(a;L, m)j2°(xL, m), where 
ripQ (xl , rn) is the number of conflgurations for the flu- 
orophore/quencher pair for a given bubble size and po- 
sition. To demonstrate this effect assume for simplic- 
ity that each bps that opens up provides one unit of 
entropy, A^fq, so that flpQ = e'"^'5f«/''s , i.e., cffcc- 
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FIG. 4: Top: Mean correlation time versus versus tempera- 
ture for the three constructs in reference . The solid curves 
are for the mean correlation time Tcorr, whereas the dashed 
curves are the longest relaxation time n. Bottom: Melting 
curves. The dashed curves are experimental results. For the 
theoretical curves we used parameters from @] for 100 mM 
NaCl concentration. 



tively we increase the statistical weight u according to 
u — > li^'^SpQ/kB ^ leading to a shift in the melting curve to- 
wards lower temperatures (as seen in experiments, com- 
pare also Fig. [4]). In reality, however, we would expect 
a more intricate m-dcpendcncc of the complexions Qfq , 
for instance, it may be so that as the first bps close the 
the tag-position open up a relatively large amount of 
fluorophore/quencher entropy is released, while further 
opening of bps contributes less. 

We stress once more that all three effects will broaden 
the relaxation time spectrum. Further control experi- 
ments will be needed to obtain more precise information 
on the effects caused by the presence of the fluorophorc 
and quencher molecules. 

The activation plot in Fig. [2] was obtained from the 
construction of the scaling plot for the blinking autocor- 
relation function by the relative shift of the individual 
curves along the logarithmic time axis. The correspond- 
ing error bars were estimated from the width of the col- 
lapsed data at the midpoint {At{xT,t) = 1/2) as 20 per 
cent of the absolute value. The real experimental error 
is likely to be higher. However, it is difficult to esti- 
mate. The activation plot indicates an Arrhcnius-type 
behaviour, that is probably due to an energetic barrier 
crossing when the bp-bond establishes. 

We point out that wc here only considered the AT9 se- 
quence from 0, and not the other to constructs A18 and 
M18. The latter two constructs have 4 or more consecu- 
tive AT-bps, and it is known that such sequences assume 
the B'-conformation rather than the usual i?-structure 
[131 for which the parameters of [l] apply. In B' DNA, 
the breathing dynamics is significantly altered [l^- Fit- 
ting our model to the A18 and MIS constructs, we found 
indeed that these sequences showed more pronounced de- 



viations from our model. 

In Fig. |4]the top panel show the mean correlation time 
Tcorr = Tf{xT,T)dT = A{xT,t)dt, scc Eqs. (ini) 
and ((24|) , for the three constructs of Ref. @ ; these con- 
structs all consist of 18 consecutive AT-bps with end- 
clamps consisting of GC-pairs. The bottom panel depicts 
the probability Poq(a;T) that the bps at xt and + 1 
are open, i.e., the probability to get a fluorescence sig- 
nal. We notice that Tcorr has pronounced maxima at 
the melting transition (the point where Peq = 1/2 in 
the bottom panel). This critical slowing down at the 
melting is indeed a characteristic signature of a phase 
transition, compare Ref. [l3|- Notice that the experi- 
mental results (dashed lines) for Pcq{xT) deviate from 
the one predicted here, indicating that the fluorophore- 
quencher pair indeed has a destabilizing effect on the 
DNA helix. Also note the different melting behaviors 
of the three construct despite identical AT and GC con- 
tents predicted here as well as by experiments; this il- 
lustrates the importance of stacking interactions. Also 
notice that there is nice agreement between our theoret- 
ical results and experiments concerning the relative or- 
dering of the melting temperatures: AT9 melts first, and 
AlS last. The horizontal line (T^ax ID) in the top panel 
represents the longest relaxation time (2M -I- l)f/7r^fc~^ 
obtained from the homopolymer model of Ref. [l3, [3 
the fimit w ^ 1, ctq ^ and c = (for M = 27, length 
of the three constructs), thus giving a scaling consistent 
with the first exit of unbiased diffusion, see Ref. [l^. 



B. Bacteriophage T7 

By master equation and stochastic simulation we in- 
vestigate the promoter sequence of the T7 phage (a bac- 
teriovirus). A promoter is a sequence (often containing 
the 4 bp long TATA motif) marking the start of a gene, to 
which RNA polymerase is recruited and where transcrip- 
tion then initiates. Previous studies [2^,[2^ based on the 
Dauxois-Peyrard-Bishop model found that the the TATA 
motif is characterized by a particularly low stability and 
therefore proneness to bubble formation, although the 
statistical relevance of those data were under discussion 
[sot . We here revisit the problem of the stability and 
dynamics of the TATA motif using the necessary full set 
of stacking interactions. The T7 promoter sequence we 
investigate is 

1 20 
I I 
5 ' -aTGACCAGTTGAAGGACTGGAAGTAATACGACTC , . 

AGTATAGGGACAATGCTTAAGGTCGCTCTCTAGGAg-3 ' ' ' 
38 41 68 

whose TATA motif is marked red [H, [H]. Fig. [H] shows 
the time series of I{t) at 37°C for the tag positions 
x-T ~ 38 in the core of TATA, and xx = 41 at the sec- 
ond GG bp after TATA. Bubble events occur much more 
frequently in TATA (the TA/AT stacking interaction is 
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FIG. 5: Top: Fluorescence time series I{t) for the T7 pro- 
moter sequence, with tag position xt = 38 (red) and xt = 41 
(blue). Middle: Waiting time (V'(''")) and fluorescence sur- 
vival time {4>{t)) densities, in units of k. The data points 
(solid lines) are results from the Gillespie algorithm (master 
equation). Bottom: Mean correlation time for A = 0. All 
results are for T = 37°C and 100 mM NaCl with DNA pa- 
rameters from @]. 

particularly weak [^). This is quantified by the density 
of waiting times iP{t) spent in the I{t) = state, whose 
characteristic time scale r' = dfrip (f ) is more than 
an order of magnitude longer than at = 41. In con- 
trast, we observe similar behavior for the density of open- 
ing times (/>(t) for xt = 38 and 41, where the characteris- 
tic time is r= dff(/)(f). The solid lines are the results 
from the master equation, see subsection IIII CI showing 
excellent agreement with the Gillespie results. Notice 
that whereas ip(t) is characterized by a single exponen- 
tial, (/)(t) show a crossover between different regimes. For 
long times both ipi''') ^^nd 4'{t) decay exponentially as 
they should for a finite DNA stretch. As shown in the 
bottom for the parameters from Q , the variation of the 
mean correlation time Tcon- = / A{xt, t)dt obtained from 
the ME is small for the entire sequence, consistent with 
the low sensitivity to the sequence of 0(t). However, 
note the even smaller variation predicted for the param- 
eters of @ , indicating that the stability parameters of [1| 
are more sequence sensitive compared to previously used 
values @. We speculate that the recurrence frequency 
of bubble events may be a clue in the understanding of 
transcription initiation: If the protein, that is supposed 
to bind to the specific site, senses a time-averaged energy 
landscape, the significantly more frequent bubble events 
at TATA may trigger its binding and thus trigger tran- 
scription initiation. 

Fig. [5] shows the equilibrium probability that the bps 
[xt — a, + A] are open, as necessary for fluorescence 
to occur. We plot data obtained from the zeroth mode 
(an ME eigenvalue problem always has one zero eigen- 
value, the corresponding eigenvalue is the equilibrium 
probability [l^l) of the ME together with the time avcr- 
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FIG. 6: Probability to have at least the bps [xt — A,xt + A] 
open (as assumed to be necessary for fluorescence to occur) 
for different tag positions xt, and for A = and 2. The 
results from master equation results and Gillespie algorithm 
show excellent agreement. Same parameters as in Fig. [5] 

age from the stochastic simulation (G), finding excellent 
agreement. Whereas for A = several segments show in- 
creased tendency to denaturation, for the case A = 2, one 
major peak is observed; the data from coincide pre- 
cisely with TATA. For comparison the equilibrium prob- 
ability obtained using DNA stability data from Q have 
their maximum peak upstream. Analysis for various A 
indicate best discrimination of the TATA sequence be- 
ing open for A = 2. Biologically, this finding is signifi- 
cant, as it corresponds to the probability for simultaneous 
opening of the whole TATA motif. For future FCS or 
energy transfer experiments investigating the relevance 
of denaturation-induced facilitation of transcription ini- 
tiation, it therefore appears important to optimize the 
A-dependence for best resolution, e.g., by adjusting the 
linker lengths of fluorophore and quencher. This could, 
in principal, be experimentally achieved as, assuming a 
circular bubble of 5 open bps with bp-bp distance 3.4A, 
the distance between fluorophore and quencher on bub- 
ble opening increases by 6-7 A, the same magnitude as 
the Forster transfer radius. 

In Fig. [5] we compare the opening probabilities to the 
values for a random sequence, for which we chose the free 
energies such that the content of AT and GC bps is 50:50. 
Then, we define 

Ehb, random = ehb,AT/2 + ehb,Gc/2, (26) 

for the hydrogen bonding, and 
1 , 

est, random — (^st,AT/TA + est,TA/AT + 2est,AT/AT 

+ est.GC/CG + est,CG/GC + 2est,GC/GC + 2est,GA/CT 
+ 2 

est,CA/GT + 2est,AG/TC + 2est,AC/TG) (27) 

for the stacking free energies. The numerical values are 
Ehb.random = 0.4 kcal/mol and est, random = -1-6 kcal/mol 
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FIG. 7: Autocorrelation function A{xT,t) for a tagged bp for 
the T7 sequence. The solid lines are master equation results 
whereas the data points represent results from the Gillespie 
scheme. Note that the autocorrelation function is sensitive 
to the DNA sequence as well as the tagging position. Same 
parameters as in Fig. [5] 
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FIG. 8: Top: Mean correlation time versus NaCl concen- 
tration for various temperatures T for the AT9 construct, 
showing a critical slowing down at the melting concentration 
(compare lower panel) . The triangles denotes the melting con- 
centration for infinitely long random AT and GC stretches. 
Bottom: Opening probability for A = 0. 



at r = 37°C and 100 mM NaCl. Inserting both into the 
expressions for the partition factor ([T]) and the rates t, 
the dashed lines in Fig. [S] are obtained. Thus, the peaks 
in the opening probabilities are distinctly significant. 

In order to illustrate the breathing dynamics of the 
T7 sequence using experimentally measurable quantities, 
Fig.[7]shows the autocorrelation functions (see subsection 
lIVBp . for four different tag positions xt (same param- 
eters as above) within the promoter region. Both the 
Gillespie approach as well as the master equation were 
used and compared; excellent agreement between them 
are found. The autocorrelation function for the tagged 
bp decays faster if positioned in a GC-rich region than 
in an AT-rich region. Comparing with Fig. [5] it should 
be possible to resolve the different decay times of the 
autocorrelation function experimentally. 

C. Nanosensing applications 

In Ref. Q the DNA Watson-Crick and stacking param- 
eters were obtained for different NaCl concentrations, al- 
lowing us to study the effect of salt concentration on the 
breathing dynamics and equilibrium properties of DNA. 
Fig. [5] shows the dependence of the mean correlation time 
Tcorr and the equilibrium opening probability Pc^{xt) for 
the AT9 sequence on salt concentration C and tempera- 
ture T, using the same tagging position as in subsection 
IIV Al We point out that the mean correlation time is 
directly accessible in experiments. Note the logarithmic 
axis. The triangles denote the melting concentration of 
infinitely long random AT and GC stretches, respectively 
(from [8|). The maxima of the Tcorr curves signify the 
critical slowing down of the autocorrelation at the phase 
transition as before; note that the maxima coincide with 



the melting concentrations in the bottom panel. The 
dashed line (xmax 2D) corresponds to the longest relax- 
ation time obtained numerically from the ME; it agrees 
well with Tcorr closc to the maximum (equivalently for 
the other T), indicating that at melting there is a sin- 
gle (slow) relevant relaxation mode. The horizontal line 
(Tmax ID) represents the analytically obtained longest 
relaxation time (2M -|- l)^/7r^fc~^ for a homopolymer 
model, compare Ref. 

The predicted variation with C and T shown in Fig. [S] 
is significant. Thus, different solvent conditions such as 
temperature or salt alter the opening probability of the 
DNA construct, and therefore the blinking activity. For 
a fixed salt concentration, for instance, a higher temper- 
ature would therefore lead to more frequent, and longer 
blinking events, such that one could measure the effect by 
both recording individual blinking events and integrating 
the blinking signals. As shown in our analysis, the sta- 
bility parameters are sufficiently sensitive to externally 
detect changes of these parameters. Note that a DNA 
construct of 30 bps roughly corresponds to a length of 
lOnm. Such nanoprobes would easily fit into nanochan- 
nels, small lipid vesicles, or microdishes in gene arrays. 
We therefore propose to investigate in more detail the 
suitability of DNA-breathing constructs as nanosensors 




V. CONCLUSIONS 

In this study we considered the bubble breathing dy- 
namics in a heteropolymer DNA-region characterized by 
statistical weights Ust{x) for disrupting a stacking inter- 
action between neighboring bps, and the weight Uhb(2;) 
for breaking a Watson-Crick hydrogen bond {x labels dif- 
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ferent bps), as well the bubble initiation parameter (the 
ring-factor) ^. For that purpose, we introduced a (2 + 1)- 
variable master equation governing the time evolution of 
the probability distribution to find a bubble of size m 
with left fork position xl at time as well as a com- 
plementary Gillespie scheme. The time averages from 
the stochastic simulation agree well with the ensemble 
properties derived from the master equation. We calcu- 
late the spectrum of relaxation times, and in particular 
the experimentally measurable autocorrelation function 
of a tagged bp is obtained. All parameters in our model 
are known from recent equilibrium measurements avail- 
able for arbitrary temperature and NaCl concentration, 
except for the rate constant k for (un)zipping that is 
the only free fit parameter. We note that the value for 
zipping rate obtained from the fluorescence correlation 
studies is significantly lower than from NMR experiments 
[3^. The difference may stem from the higher temper- 
atures and longer AT sequences probed in the fluores- 
cence experiments. However, a perturbing effect of the 
fluorophore-quencher pair in the FCS approach cannot 
be excluded. For a better understanding of fc, a more de- 
tailed microscopic modeling and additional experimental 
study are needed. 

We applied recent DNA stability data from based 
on separation of hydrogen bond and stacking energies. A 
distinct feature of these parameters is the low stacking 
in an TA/AT pair of bps, translating into a pronounced 
instability of the TATA motif, as shown for the T7 pro- 
moter sequence. We demonstrated that the probability 
of simultaneous opening of a stretch of the size of 4 to 
5 bps well discriminates the TATA motif from the other 
positions along the promoter sequence, reflecting its bi- 
ological relevance. This demonstrates that single DNA 
fluorescence spectroscopy experiments can likely be used 
to investigate in more details the role of the interplay 
between TATA-breathing. TATA-box binding proteins, 
and transcription initiation. Regarding the mechanism 
how TATA may guide this initiation we speculate that it 
is not primarily the bubble lifetime (much shorter than 



the timescale of typical conformational changes of pro- 
teins) but the recurrence frequency of bubble events that 
triggers the protein binding. 

We note that there exists also a Langevin equation ap- 
proach to DNA-breathing, the Dauxois-Peyrard-Bishop 
model [3^ . [ssj . with seven free parameters. Values of 
these parameters were assigned by comparison to exper- 
imental melting curves for three different short DNA se- 
quences obtained for rather specific solvent conditions in 
[3q . In particular, stacking interactions were taken to 
be independent of bp sequence [36| . In view of the direct 
measurement of the stacking free energy in Q under var- 
ious conditions, it would be desirable to modify the DPB 
model to accommodate for the full set of new stability 
parameters. 

We expect this study to encourage furthcrgoing in- 
vestigations on the theoretical understanding of DNA- 
breathing and the experimental possibilities to obtain 
detailed sequence and stability information of DNA 
and its interactions with binding proteins from DNA- 
breathing dynamics. We furthermore point out the pos- 
sibility to use the results of this study for designing a 
small fiuorophore/quencher-dressed DNA construct for 
nanosensing applications in nanochannels, vesicles or mi- 
crodishes. 
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